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ABSTRACT 

The  seismic  discrimination  between  weapons  tests,  industrial  explosions,  natural  earthquakes,  and 
other  types  of  seismic  events  is  limited  by  an  incomplete  understanding  of  the  basic  physics  of  the  source 
processes.  The  problem  of  quantifying  how  the  generation  of  elastic  waves  is  related  to  such  factors  as 
mode  of  energy  release,  source  dimensions,  material  properties,  and  stress  concentrations  is  still  unsolved 
for  the  most  part.  While  a  complete  solution  of  this  complicated  problem  may  not  be  possible  at  the 
present  time,  an  examination  of  the  similarities  and  differences  between  various  types  of  seismic  sources 
could  prove  useful  in  the  discrimination  task.  This  was  one  of  the  objectives  of  this  research  effort. 

A  starting  point  for  the  comparison  of  the  different  types  of  seismic  sources  is  the  large  amounts  of 
empirical  data  that  are  summarized  in  various  seismic  scaling  relationships,  such  as  moment  versus  yield 
for  contained  explosions  and  stress  drop  versus  moment  for  natural  earthquakes.  By  converting  these 
scaling  relationships  for  different  types  of  sources  to  a  common  form,  it  is  possible  to  assess  whether 
there  are  significant  differences  in  the  basic  source  processes.  For  instance,  both  contained  explosions 
and  natural  earthquakes  are  consistent  with  a  linear  relationship  between  seismic  scalar  moment  and 
source  energy,  with  the  proportionality  factor  being  identical  if  the  efhciency  of  the  earthquake  source  is 
taken  to  be  about  0.15%.  A  recent  study  of  the  energy  required  to  generate  pseudotachylites  suggests 
that  this  value  of  efhciency  is  not  unreasonable.  Another  approach  to  the  comparison  of  seismic  sources 
is  the  consideration  of  energy  density  in  the  source  region,  which  is  related  to  the  total  energy  of  the 
source  and  its  spatial  dimensions.  Here  again  recently  acquired  data  suggest  that,  provided  a  distinction 
is  made  between  the  initial  source  dimension  and  the  effective  source  dimension,  the  energy  densities  of 
different  types  of  sources  are  more  similar  than  previously  thought.  This  leads  to  a  consideration  of  such 
matters  as  the  interaction  between  the  seismic  source  and  the  material  properties,  the  degree  to  which  the 
material  properties  are  modihed  by  the  source,  and  the  possibility  of  secondary  sources.  These  aspects 
of  the  problem  have  been  examined  by  simulating  the  stress  field  near  different  types  of  seismic  sources, 
comparing  the  results,  and  exploring  the  effect  of  various  source  parameters. 
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OBJECTIVES 

The  scope  of  this  research  effort  is  concerned  with  the  development  of  an  improved  understanding  of 
methods  for  locating  and  characterizing  seismic  events  in  a  heterogeneous  earth.  The  general  objective  is 
to  investigate  how  improved  models  for  the  generation  and  propagation  of  elastic  waves  can  help  in  the 
evaluation  of  methods  and  models  currently  being  used.  The  work  statement  for  the  research  contains 
two  primary  tasks; 

1.  Continue  the  development  of  analytical  and  numerical  techniques  capable  of  modeling  the  physical 
processes  that  cause  the  generation  of  elastic  waves  by  seismic  sources.  Also  continue  the  develop¬ 
ment  of  three-dimensional  models  of  the  earth’s  velocity  structure  which  have  sufhcient  resolution 
to  provide  accurate  locations  of  seismic  events  at  regional  and  teleseismic  distances. 
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2.  The  modeling  capability  will  be  validated  against  existing  data  bases  and  any  new  data  which 
can  be  acquired  by  taking  advantage  of  targets  of  opportunity  to  perform  broad  band  recording 
experiments. 

This  paper  reports  on  work  that  has  bee  done  during  the  past  three  years  on  scaling  relationships 
for  seismic  sources,  with  particular  emphasis  on  such  relationships  for  small  earthquakes.  Developing 
a  more  complete  understanding  of  these  scaling  relationships  for  both  contained  explosions  and  natural 
earthquakes  has  the  potential  to  provide  a  more  comprehensive  theoretical  foundation  for  the  methods 
that  are  currently  being  used  to  monitor  a  Comprehensive  Nuclear- Test-Ban  Treaty  and  perhaps  uncover 
new  methods. 

RESEARCH  ACCOMPLISHED 

Haskell  (1964)  introduced  a  model  of  an  earthquake  as  a  rupture  propagating  over  a  finite  section 
of  a  fault  and  provided  convenient  analytical  formulas  for  the  radiated  elastic  waves.  A  characteristic 
of  this  model  in  its  most  simple  form  is  a  very  uniform  release  of  displacement  and  stress  on  the  fault 
surface.  Haskell  recognized  this  as  a  deficiency  and  suggested  that  some  form  of  heterogeneity  was  needed 
in  the  distribution  of  displacement  in  either  time  or  space.  Haskell  (1966)  and  Aki  (1967)  attempted  to 
model  this  heterogeneity  by  describing  the  acceleration  and  velocity,  respectively,  of  fault  displacement 
in  terms  of  random  functions  that  were  correlated  over  only  limited  portions  of  the  fault  plane.  Brune 
(1970)  considered  the  time  history  of  displacement  on  the  fault  and  provided  spectral  models  relating 
the  frequency  content  of  radiated  elastic  waves  to  dimensions  and  stress  present  on  the  fault.  These 
models  for  an  earthquake  formed  the  basis  for  the  interpretation  of  large  amounts  of  observational  data 
in  terms  of  earthquake  source  parameters,  but  the  emphasis  was  primarily  upon  average  properties  and 
not  much  attention  was  given  to  the  heterogeneity  originally  explored  by  Haskell  and  Aki.  Kanamori  and 
Anderson  (1975)  summarized  many  of  these  empirical  data  and  used  them  to  justify  scaling  laws  that 
related  such  parameters  as  fault  dimension,  stress  drop,  moment,  magnitude,  and  energy.  Emerging  from 
studies  of  this  type  was  a  generally  held  conviction  that  the  stress  drops  of  most  shallow  earthquakes  was 
independent  of  size  and  generally  fell  in  the  range  of  10  to  100  bars.  Implicit  in  most  of  these  studies  was 
the  assumption  that  the  fault  was  relatively  uniform  in  terms  of  material  properties  and  energy  release. 

A  rather  different  approach  to  the  study  of  heterogeneity  on  the  fault  surface  has  emerged  rather 
unexpectedly  from  studies  of  small  repeating  earthquakes  along  the  San  Andreas  fault  system  in  Califor¬ 
nia.  Nadeau  and  Johnson  (1998)  showed  how  measurements  of  moment  release  rates  of  these  repeating 
earthquakes  could  be  used  to  construct  scaling  relationships  between  seismic  moment,  repeat  time,  fault 
slip  and  fault  dimension.  Their  method  of  interpreting  the  data  does  not  require  the  use  of  any  model  for 
the  earthquake  process  and  uses  the  single  assumption  that  slip  on  the  fault  at  depth  is  closely  related 
to  the  measured  slip  at  the  surface.  The  study  of  Nadeau  and  McEvilly  (1999)  provides  strong  support 
for  this  assumption  by  showing  that  second  order  variations  of  slip  in  both  time  and  space  are  strongly 
correlated  between  the  surface  and  the  depth  of  the  seismogenic  zone.  The  basic  scaling  relationship 
obtained  by  Nadeau  and  Johnson  (1998)  relates  repeat  time  of  earthquakes  T  to  scalar  seismic  moment 
according  to  the  formula 

T  oc  MJ  (1) 

where  the  fraction  ^  is  a  close  approximation  to  the  numerical  value  of  0.17  which  was  obtained  by  regres¬ 
sion  analysis.  This  result  can  also  be  interpreted  in  terms  of  scaling  relationships  for  fault  displacement 
u  and  fault  area  A  of  the  form 

u  oc  Me?  (2) 

and 

A  oc  m!  (3) 

When  these  relationships  are  interpreted  in  terms  of  stress  drops  they  lead  to  values  in  the  GPa  range, 
which  is  not  impossible  for  crustal  rocks  (Sammis  et  ah,  1999),  but  much  larger  than  expected  values  in 


the  1  to  10  MPa  range.  This  result,  together  with  the  fact  that  the  repeating  earthquakes  occupy  a  very 
small  fraction  of  the  fault  and  are  surrounded  by  regions  which  appear  to  be  creeping  and  thus  quite  weak, 
argues  for  a  fault  which  is  highly  heterogeneous  in  its  strength  properties.  Alternative  interpretations 
have  been  recently  presented  (Anooshehpoor  and  Prune,  1998,  2000;  Sammis  and  Rice,  2000;  Beeler, 
2000)  which  attempt  to  avoid  the  need  for  high  stress  drops  and  thus  do  not  require  the  heterogeneous 
fault  surface. 

While  one  of  the  strengths  of  the  study  by  Nadeau  and  Johnson  (1998)  was  that  it  did  not  use  a 
source  model  that  presumed  either  a  homogeneous  or  heterogeneous  fault  surface,  a  full  interpretation  of 
the  results  does  require  that  some  type  of  model  be  used.  Here  we  introduce  such  a  model  and  use  this 
model  to  interpret  the  data  in  terms  of  physical  properties  and  processes  on  the  fault  surface.  The  model 
introduced  here  begins  with  the  assumption  that  strength  on  the  fault  surface  is  highly  heterogeneous  and 
attempts  to  show  that  such  a  model  is  consistent  with  the  scaling  relationships  presented  in  Nadeau  and 
Johnson  (1998).  Only  the  quasi-static  part  of  the  earthquake  process  is  considered  here,  which  extends 
up  to  the  time  that  rupture  begins. 

Consider  a  planar  fault  surface  lying  in  the  xy  plane  and  let  the  slip  across  this  surface  be  denoted  by 

d(a;,  y)  =  d^{x,  y)x  -f  dy{x,  y)y  (4) 

Without  loss  of  generality  we  will  assume  that  the  motion  on  the  fault  is  in  response  to  a  loading  traction 
Tx  that  acts  only  in  the  x  direction.  Let  the  slip  that  would  occur  if  the  fault  were  uniformly  weak  and 
creeping  uniformly  be  denoted  by  d^.  Then  the  displacement  deficit  on  the  fault  can  be  defined  as 

u(a;,  y)  =  dc  -  d(a;,  y)  (5) 


At  points  where  the  fault  slip  is  keeping  up  with  the  applied  traction  the  displacement  deficit  will  be 
zero. 

The  fault  is  considered  to  be  sufhciently  weak  so  that  it  continually  creeps  as  tectonic  stress  is  applied 
except  on  isolated  patches  which  we  will  call  asperities.  The  asperities  are  small  areas  of  the  fault  which 
are  much  stronger  than  the  surrounding  regions  and  able  to  withstand  the  tectonic  stress  until  a  stress 
limit  is  reached  and  rupture  occurs.  These  asperities  will  be  assumed  to  be  circular  with  radius  and 
area  =  Trr^.  In  the  static  problem  being  considered  here  we  only  consider  the  situation  up  to  the 
instant  just  before  the  asperities  rupture,  so  there  will  be  no  movement  on  these  asperities.  Thus  the 
displacement  deficit  u  will  be  a  maximum  at  the  asperities. 

Focus  now  upon  a  region  of  the  fault  where  a  group  of  asperities  are  close  enough  together  so  that 
they  strongly  interact.  Assume  that  this  group  of  asperities  is  roughly  in  the  form  of  a  circle  with  radius 
To  and  area  Aa  =  Trr^.  The  displacement  deficit  will  be  the  same  on  all  asperities  and  denoted  by  Ua- 
The  asperities  are  assumed  to  be  close  enough  together  so  that  they  prevent  almost  all  motion  on  the 
intervening  areas  between  the  asperities,  and  thus  the  displacement  deficit  is  very  close  to  Ua  over  the 
entire  area  Aa-  A  solution  for  elastostatic  problem  just  described  can  be  found  in  Westmann  (1968).  The 
stress  on  the  region  Aa  where  the  displacement  is  Ua  is  given  by 


CTxzir) 


AflUa  1 

(2  —  /z)7r  _  J.2 


(6) 


where  =  x^  +  y"^ ,  fi  is  the  shear  modulus,  and  /z  is  Poisson’s  ratio.  The  stress  on  the  weak  part  of  the 
fault  outside  Aa  is  assumed  to  be  zero  in  this  solution.  However,  the  displacement  caused  by  the  group 
of  asperities  does  extend  outside  Aa  and  in  the  region  r  >  ra  is  given  by 


u^{r,  e) 
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where  9  is  the  polar  angle  in  the  xy  plane.  We  see  that  the  group  of  asperities  is  surrounded  on  the  fault 
plane  by  a  displacement  shadow  that  decays  as  r~^  at  distances  much  greater  than  r^.  Note  that  the  first 


term  in  the  expression  for  is  about  7  times  larger  than  the  second  term  in  Mj;  and  the  only  term  in 
Uy.  Thus  it  will  be  sufficient  in  what  follows  to  assume  that  the  displacement  deficit  within  the  shadow 
is  given  approximately  by 


Ux{r) 


UaTa 


r 


(9) 


Next  we  will  assume  that  the  displacement  shadow  of  the  group  of  asperities  is  terminated  when 
the  displacement  falls  below  a  critical  level  denoted  by  u^-  One  can  appeal  to  state  variable  friction 
theories  to  justify  the  existence  of  such  a  critical  level  of  displacement.  The  radius  of  the  displacement 
shadow  Tc  is  then  given  by 

rc  =  —Ta  (10) 


The  geometrical  picture  of  the  situation  on  the  fault  plane  is  now  complete,  and  it  is  shown  in  Figure  1. 
It  will  be  assumed  that  Ua  Mc  and  it  follows  that  Ta- 

Were  all  of  the  asperities  to  rupture  and  all  of  the  displacement  deficit  released,  the  scalar  moment 
of  such  an  event  would  be  approximately 


Mo  = 
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From  the  expression  for  stress  on  the  asperity  patch  (6)  it  is  clear  that  this  stress  in  not  uniform 
and  is  greatest  near  the  edge  of  the  patch.  Thus,  assuming  the  individual  asperities  are  similar  in  their 
properties,  one  would  expect  the  failure  of  the  asperity  patch  to  begin  with  the  failure  of  one  of  those 
closest  to  border  r^.  This  is  consistent  with  the  results  of  Das  and  Kostrov  (1983,  1986)  who  found  that 
failure  of  an  asperity  began  on  its  border  and  then  progressed  toward  the  center.  The  average  stress  on 
one  of  these  border  asperities  can  be  evaluated  to  yield 

(12) 


Ja, 


y)dxdy  =  C 


where 


32/r 

377^(2  —  v) 


(13) 


and  it  has  been  assumed  that  Va  To  in  evaluating  the  integral.  Note  that  this  result  is  similar  to  that 
of  Sammis  and  Rice  (2000),  who  used  it  to  argue  for  a  very  weak  asperity.  Here  we  use  it  to  argue  that 
when  a  group  of  asperities  are  combined  together  in  close  proximity  on  a  fault,  they  interact  in  such  way 
as  to  reduce  the  stress  felt  by  any  individual  asperity.  In  effect,  the  displacement  field  is  smoothed  by 
the  interaction  of  the  asperities,  thus  reducing  displacement  gradients,  strain,  and  stress.  Assuming  that 
dmax  is  a  constant  that  controls  the  initiation  of  rupture  on  the  asperity  patch,  we  have 


^max  , - 
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and  then 
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Figure  1:  Schematic  sketch  of  the  earthquake  model  consisting  of  a  patch  of  interacting  strong  asperities. 
The  upper  panel  shows  a  group  of  asperities  having  radius  with  a  strong  interacting  region  of  radius 
Tfl  and  surrounded  by  a  displacement  shadow  that  extends  out  to  a  radius  r^.  The  lower  panel  shows  the 
disiplacement  deficit  as  a  function  of  radial  distance  from  the  center  of  the  asperity  patch. 


The  next  step  is  to  consider  how  the  results  given  above  change  when  the  size  of  the  asperity  patch 
changes.  For  this  it  is  convenient  to  introduce  the  number  n,  which  is  taken  to  be  the  number  of  small 
asperities  in  the  asperity  patch.  If  the  asperities  were  to  add  linearly  with  no  interaction,  then  one  would 
expect  that 

Aa  =  Trrln  (16) 

which  would  imply  that 

ra{n)  =  roTii  (17) 

If  the  asperities  were  to  add  with  some  interaction  between  them,  then  we  might  expect  a  more  general 
relationship  of  the  form 


ra{n)  =  roTi'"  (18) 

where  k>  but  is  otherwise  undetermined  at  this  point. 

Now  consider  what  happens  when  the  number  of  asperities  in  the  patch  is  increased  by  one.  Then  we 
have 

ra{n  +  !)  =  («  +  l)^ro  «  ra{n){l  +  ^)  (19) 


r,(n  +  1)  «  r,(n)(l  +  |-) 

Ua{n  +  1)  fS  Ua{n){l  +  ^) 
and  finally 

Mo{n  +  l)^Mo{n){l  +  —) 

n 

From  this  last  equation  we  can  form  the  differential  equation 

^Mo{n)  =  Mo{n)  — 
an  n 

which  has  the  solution 
M„{n)  = 

The  same  process  produces 

ra{n)  =  ra{i)n'^  =  rou'^ 

z'  \  /'ix  -K  ^max  2  -k 

rc{n)  = 


(20) 

(21) 

(22) 

(23) 

(24) 

(25) 

(26) 


Ua{n)  =  Ua{l)n^  = 
and  from  these  it  is  possible  to  obtain 
M„{n)  =  27r/r^r4n3- 

O  Uq 


(27) 


(28) 


These  last  four  equations  form  a  consistent  set  of  parametric  equations  for  all  of  the  geometrical 
properties  of  this  asperity  model  of  the  earthquake  process.  They  are  all  expressed  in  terms  of  the 
parameter  n,  which  was  introduced  as  the  number  of  asperities  involved,  but  may  be  considered  more 
generally  as  a  surrogate  for  the  size  of  the  earthquake. 


Using  the  expressions  listed  above  for  ra(n)  or  rc(n)  one  can  calculate  the  areas  of  the  asperity  patch 
or  the  displacement  shadow,  respectively.  However,  neither  of  these  corresponds  to  the  area  that  was 
calculated  by  Nadeau  and  Johnson  (1998).  If  we  define 

M„{n)  =  fxA(n)ua(n)  (29) 

then  A(n)  is  an  effective  area  that  corresponds  to  the  area  calculated  by  Nadeau  and  Johnson  (1998).  It 
follows  that 


A{n)  =  A„(l)n^'^  =  (30) 

O  Uq 

These  results  for  displacement,  area,  and  scalar  moment  can  be  expressed  in  terms  of  the  scaling  rela¬ 
tionships 

Ua  OC  (31) 

A  OC  M®/®  (32) 


Assuming  that  the  displacement  on  the  asperities  Ua  accumulates  at  a  steady  rate  in  time  t  then  we  also 
have 


t  OC  M^® 


(33) 


This  is  the  basic  observational  result  presented  by  Nadeau  and  Johnson  (1998).  Note  that  all  of  these 
scaling  relationships  are  independent  of  the  parameter  n  and  the  parameter  k,  which  together  control  the 
size  of  the  area  where  there  is  strong  interaction  between  the  asperities. 

Stress  drop  is  a  parameter  that  is  commonly  used  to  characterize  earthquakes,  with  conventional 
thought  being  that  stress  drop  is  independent  of  earthquake  size  and  in  the  range  of  10  to  100  bars 
(Kanamori  and  Anderson,  1975).  In  the  asperity  model  presented  in  this  paper  stress  is  highly  heteroge¬ 
neous,  ranging  from  infinity  (in  a  mathematical  sense)  on  the  border  of  a  single  asperity  to  a  very  small 
value  in  the  displacement  shadow.  This  makes  it  difhcult  to  characterize  the  stress  with  a  single  number, 
so  several  available  options  will  be  discussed.  The  only  stress  that  is  critical  to  the  theory  is  the  average 
stress  that  causes  an  asperity  near  the  edge  of  the  patch  to  fail  (Tmax,  which  was  assumed  to  be  a  constant 
material  property.  This  can  be  compared  with  the  average  stress  felt  by  a  single  isolated  asperity 
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7r(2  —  ly) 
From  this  we  have 
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m.  a ; 
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which  shows  clearly  how  interaction  of  a  group  of  asperities  can  reduce  the  maximum  stress.  Note  that 
it  was  assumed  that  n  1  in  deriving  ^max  so  this  result  does  not  give  the  correct  answer  as  n  ^  1. 
Stress  drop  is  commonly  calculated  using  the  formula  (Kanamori  and  Anderson,  1975) 


Act  = 


l-Kfl  u 
16  r 


(36) 


where  u  is  the  average  slip  over  a  circular  crack  of  radius  r.  Using  Ua{n)  for  u  and  ra{n)  for  r  yields 


AcTfl 
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Still  another  possibility  is  to  use  the  formula  (Kanamori  and  Anderson,  1975) 


Act  = 


16 


(38) 
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Using  rc(n)  for  r  yields 


Acr^  = 


7  Mo{n)  _  28 


16  rc(n)3  37r(2  -  i/)  (Tmaxri 
The  stress  drops  calculated  in  Nadeau  and  Johnson  (1998)  correspond  to  taking  Ua{n)  for  u  and 


A{n)/Tr  for  r  and  this  yields 


Ttt//  Ua{n)  77r^\/3//(2  —  v) 


16 


128 


/  A{n)/Tr 

If  this  result  is  scaled  against  moment  we  obtain 
Ad"  oc 

which  agrees  with  the  empirical  result  Nadeau  and  Johnson  (1998)  obtained. 
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CONCLUSIONS  AND  RECOMMENDATIONS 

The  model  of  an  earthquake  that  has  been  presented  here  is  based  on  the  concept  of  small  strong 
asperities  that  resist  motion  on  the  fault.  When  a  group  of  these  asperities  are  close  enough  together  so 
that  they  can  strongly  interact,  an  asperity  patch  is  formed  which  is  stronger  than  any  of  the  individual 
asperities  and  so  a  greater  amount  of  tectonic  displacement  can  accumulate  before  failure  occurs.  The 
asperity  patch  is  surrounded  by  a  displacement  shadow  where  creep  displacement  on  the  weak  part  of 
the  fault  is  retarded  because  of  the  presence  of  the  asperity  patch.  The  asperity  patch  ruptures  and  an 
earthquake  occurs  when  the  stress  on  one  of  the  asperities,  most  likely  near  the  border  of  the  patch,  reaches 
its  strength  limit  and  fails,  which  causes  a  cascade  failure  of  all  the  other  asperities  in  the  patch.  At 
this  time  the  entire  displacement  deficit  associated  with  the  asperity  patch,  that  within  both  the  asperity 
patch  and  the  surrounding  displacement  shadow,  is  released  and  contributes  to  the  scalar  moment  of  the 
earthquake.  This  model  depends  upon  only  a  few  parameters,  the  radius  of  the  individual  asperities  Tq, 
the  critical  level  of  displacement  that  determines  the  limit  of  the  displacement  shadow  u^,  the  stress  level 
at  which  an  individual  asperity  will  fail  (Tmax,  the  number  of  asperities  in  an  asperity  patch  n,  and  the 
exponential  factor  that  controls  how  the  radius  of  the  asperity  patch  grows  with  the  number  of  asperities 
K.  What  is  rather  remakable  is  that  the  scaling  relationships  for  displacement,  area,  and  time  to  failure 
as  a  function  of  moment  do  not  depend  upon  any  of  these  parameters.  They  are  direct  consequences  of 
analytic  elastostatic  solutions  for  the  exterior  crack  problem.  Another  characteristic  of  this  model  is  that 
it  involves  a  very  heterogeneous  stress  field  that  is  difhcult  to  characterize  in  a  meaningful  way  with  a 
single  parameter  such  as  stress  drop. 
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